function [h, K] = mimpresp (N, D, t)
%
% h = mimpresp (N, D, t)  
%
% Calculate impulse response for transfer function with numerator
% N and denominator D for times t.  An error results if the impulse
% response contains impulses and/or derivatives of impulses.
%
% Note this is numerically ill-conditioned if there are closely
% spaced poles.
%
% See also mresponse.  
  
% Michael Hayes, UCECE

  if 0  
    ss = tf (N, D);
    % This aborts if transfer function not strictly proper and
    % quietly gives the wrong answer when numerator degree is
    % the same as denominator degree.
    h = impulse (ss, t);
    return
  end
  
  % Remove common factors of s.  This is pole-zero cancellation.
  while N(end) == 0 & D(end) == 0
    N = N(1:end - 1)
    D = D(1:end - 1)
  end
  
  [R, P, K] = residue (N, D);
  
  if ~isempty (K)
    error (['Impulse response has impulses and/or derivatives of ' ...
            'impulses.'])
  end

  % Sort the poles.
  [PS, I] = sort (P);
  RS = R(I);

  NP = numel (P);
  
  h = zeros (size (t));
  y = zeros (size (t));  

  dt = t(2) - t(1);

  n = 1;
  while n <= NP
    p = PS(n);
    % Find out many times this pole is repeated.
    mp = find (PS == p);
    NM = numel (mp);
    
    for m = 1 : NM
      h = h + RS(n) / factorial (m -1) * t.^(m - 1) ...
          .* exp (PS(n) * t) .* (t >= 0);        
      n = n + 1;
    end
  end

  h = real (h);
  

